enrichMiR benchmark results

suppressPackageStartupMessages({
  library(ggplot2)
  library(plgINS)
  library(DT)
  library(viridis)
  library(RColorBrewer)
  library(data.table)
  library(dplyr)
  library(cowplot)
  library(matrixStats)
  library(RUVSeq)
})

Allocation

# input file: Bartel DEA SE
input <- "/mnt/schratt/enrichMiR_benchmark/results1/aggregated.benchmark.rds"

Loading

# import
df.all <- readRDS(input)

test <- split(df.all, df.all$dataset)
df <- test$bartel.hek
dea.names <- as.character(unique(df$treatment))
props.all <- levels(df$prop)

Plots

# plots
for(i in c("detPPV","FP.atFDR05","log10QDiff","TP.atFDR05")){
  cat("Score analysis: ", i,"\n\n")
  ## boxplot
  # print( 
  #   ggplot(df, aes(x=prop, y=df[[i]], fill=method)) + geom_boxplot() + ylab(i) 
  #   )
  if(i!="TP.atFDR05"){
    ## density plot: score distribution per permutation
    print( 
      ggplot(df, aes(log10(df[[i]]))) + geom_density(aes(fill=method), alpha=0.9)  + facet_wrap(~prop) + xlab(i)
      )
    ## violin plot
    if(!i %in% c("FP.atFDR05","log10QDiff")){
      print( 
        ggplot(df, aes(x=prop, y=df[[i]], fill=method))  + geom_violin() + facet_wrap(~method) +
          guides(fill=FALSE) + ylab(i) + geom_jitter(shape=16, position=position_jitter(0.2))
      )
    }
  }
  ## boxplot FP
  if(i=="FP.atFDR05"){
    print(
      ggplot(df, aes(x=method, y=df[[i]], fill=method)) + geom_boxplot() +
            guides(fill=FALSE) + ylab(i) + coord_flip() + scale_x_discrete(limits = rev(levels(df$method)))
          )
  }
  ## violin plot median
 print(
   ggplot(df, aes(x=prop, y=df[[i]], fill=method))  + geom_violin(draw_quantiles = 0.5) + facet_wrap(~method) +
    guides(fill=FALSE) + ylab(i)
  ) 
  ## mean curves plot
  if(i!="FP.atFDR05"){
    print(
      ggplot(df, aes(x=as.integer(df$prop), y=df[[i]], color=method)) + geom_smooth(alpha=.1) + scale_x_continuous(labels=unique(props.all)) + 
        ylab(i) + xlab("prop")
      )
  }
}

## Score analysis:  detPPV

## Score analysis:  FP.atFDR05

## Score analysis:  log10QDiff

## Score analysis:  TP.atFDR05

# FP-TP plot (at FDR .05)
## get average number of FP & TP at FDR .05 for each combination of permutation & method
df.agg <- aggregate(df[,c("FP.atFDR05","TP.atFDR05")], by=df[,c("prop","method")], FUN=mean)
## plot
ggplot(df.agg, aes(x=FP.atFDR05, y=TP.atFDR05, color=method, shape=prop, group=method)) + geom_line() + geom_point(size=3) + 
  geom_label(data=df.agg %>% filter(prop==50), aes(label=method), position=position_jitter(width=3,height=.03))

scores over the datasets

# mean scores over all datasets; this shows which datasets were difficult to handle for the tests
ggplot(df, aes(x=prop, y=detPPV, fill=treatment))  + geom_violin(draw_quantiles = 0.5) + facet_wrap(~treatment) +
      guides(fill=FALSE)

mean.datasets <- lapply( unique(df$treatment), FUN=function(x) 
  sapply( unique(df$prop), function(y) 
    round( mean(df$detPPV[df$treatment==x & df$prop==y], na.rm=TRUE), 3) 
    )
  )
names(mean.datasets) <- dea.names
for(i in dea.names){
  names(mean.datasets[[i]]) <- unique(df$prop)
}
mean.datasets

## $miR.122
## original       20       35       50 
##    0.732    0.732    0.732    0.732 
## 
## $miR.133
## original       20       35       50 
##    0.733    0.733    0.733    0.733 
## 
## $miR.138
## original       20       35       50 
##    0.694    0.730    0.700    0.696 
## 
## $miR.145
## original       20       35       50 
##    0.686    0.686    0.686    0.650 
## 
## $miR.184
## original       20       35       50 
##    0.600    0.598    0.600    0.598 
## 
## $miR.190a
## original       20       35       50 
##    0.732    0.732    0.702    0.732 
## 
## $miR.200b
## original       20       35       50 
##    0.752    0.745    0.752    0.739 
## 
## $miR.216a
## original       20       35       50 
##    0.628    0.633    0.612    0.561 
## 
## $miR.217
## original       20       35       50 
##    0.698    0.699    0.696    0.698 
## 
## $miR.219a
## original       20       35       50 
##    0.732    0.732    0.732    0.732 
## 
## $miR.375
## original       20       35       50 
##    0.707    0.707    0.707    0.707 
## 
## $miR.451a
## original       20       35       50 
##    0.687    0.687    0.672    0.709

scores over methods

# mean scores over all methods
ggplot(df, aes(x=prop, y=detPPV, fill=method))  + geom_violin(draw_quantiles = 0.5) + facet_wrap(~method) +
      guides(fill=FALSE)

mean.methods <- lapply( unique(df$method), FUN=function(x) 
  sapply( unique(df$prop), function(y) 
    round( mean(df$detPPV[df$method==x & df$prop==y], na.rm=TRUE), 3)
    )
  )
names(mean.methods) <- unique(df$method)
for(i in unique(df$method)){
  names(mean.methods[[i]]) <- unique(df$prop)
}
mean.methods

## $aREAmir
## original       20       35       50 
##    0.958    0.958    0.958    0.958 
## 
## $aREAmir2
## original       20       35       50 
##    0.958    0.958    0.958    0.958 
## 
## $EN.up
## original       20       35       50 
##    0.006    0.006    0.006    0.006 
## 
## $EN.down
## original       20       35       50 
##        1        1        1        1 
## 
## $wEN.up
## original       20       35       50 
##    0.005    0.005    0.005    0.005 
## 
## $wEN.down
## original       20       35       50 
##        1        1        1        1 
## 
## $michael.up
## original       20       35       50 
##    0.005    0.005    0.005    0.005 
## 
## $michael.down
## original       20       35       50 
##        1        1        1        1 
## 
## $MW
## original       20       35       50 
##    0.876    0.876    0.876    0.918 
## 
## $KS
## original       20       35       50 
##    0.753    0.753    0.726    0.753 
## 
## $KS2
## original       20       35       50 
##    1.000    1.000    1.000    0.958 
## 
## $GSEA
## original       20       35       50 
##    0.584    0.662    0.551    0.527 
## 
## $modSites
## original       20       35       50 
##        1        1        1        1 
## 
## $modScore
## original       20       35       50 
##        1        1        1        1 
## 
## $regmir.down
## original       20       35       50 
##    0.201    0.179    0.198    0.169 
## 
## $regmir.up
## original       20       35       50 
##    0.009    0.009    0.009    0.009 
## 
## $regmirb.down
## original       20       35       50 
##        1        1        1        1 
## 
## $regmirb.up
## original       20       35       50 
##     0.01     0.01     0.01     0.01 
## 
## $combFish.1
## original       20       35       50 
##        1        1        1        1 
## 
## $combGeom.1
## original       20       35       50 
##        1        1        1        1 
## 
## $combFish.2
## original       20       35       50 
##    1.000    1.000    0.958    0.958 
## 
## $combGeom.2
## original       20       35       50 
##    1.000    1.000    1.000    0.958

TP ratio at FDR 0.05 over methods

# sensitivity score (ratio of TP at FDR 0.05)
sens.methods <- lapply( unique(df$method), FUN=function(x) 
  sapply( unique(df$prop), function(y)
    round( sum(df$TP.atFDR05==1 & df$method==x & df$prop==y) / sum(df$method==x & df$prop==y), 3)
    )
  )
names(sens.methods) <- unique(df$method)
for(i in unique(df$method)){
  names(sens.methods[[i]]) <- unique(df$prop)
}
sens.methods

## $aREAmir
## original       20       35       50 
##    0.917    0.917    0.917    0.917 
## 
## $aREAmir2
## original       20       35       50 
##    0.917    0.917    0.917    0.917 
## 
## $EN.up
## original       20       35       50 
##        0        0        0        0 
## 
## $EN.down
## original       20       35       50 
##        1        1        1        1 
## 
## $wEN.up
## original       20       35       50 
##        0        0        0        0 
## 
## $wEN.down
## original       20       35       50 
##        1        1        1        1 
## 
## $michael.up
## original       20       35       50 
##        0        0        0        0 
## 
## $michael.down
## original       20       35       50 
##        1        1        1        1 
## 
## $MW
## original       20       35       50 
##        1        1        1        1 
## 
## $KS
## original       20       35       50 
##        1        1        1        1 
## 
## $KS2
## original       20       35       50 
##    0.917    0.917    0.917    1.000 
## 
## $GSEA
## original       20       35       50 
##    0.250    0.250    0.333    0.250 
## 
## $modSites
## original       20       35       50 
##        1        1        1        1 
## 
## $modScore
## original       20       35       50 
##        1        1        1        1 
## 
## $regmir.down
## original       20       35       50 
##        0        0        0        0 
## 
## $regmir.up
## original       20       35       50 
##        0        0        0        0 
## 
## $regmirb.down
## original       20       35       50 
##    0.833    0.833    0.833    0.833 
## 
## $regmirb.up
## original       20       35       50 
##        0        0        0        0 
## 
## $combFish.1
## original       20       35       50 
##        1        1        1        1 
## 
## $combGeom.1
## original       20       35       50 
##        1        1        1        1 
## 
## $combFish.2
## original       20       35       50 
##        1        1        1        1 
## 
## $combGeom.2
## original       20       35       50 
##        1        1        1        1